Виконавці: Данило Кравчук (завантаження даних, DESeq, Функц. аналіз) Лев Кравчук (PCA, Функц. аналіз)

GitHub

Мета

Для аналізу диференційної експресії генів було обрано клітини лінії А549 як біологічний об’єкт. Клітини лінії A549 — це клітини недрібноклітинної (NSCLC) аденокарциноми легень.

Завантаження даних

Дані було отримано з бази даних GEO (Gene Expression Omnibus). Завантажені дані описують експеримент з секвенуванням РНК (RNA-sequencing) на опромінених рентгенівським випромінюванням і контрольних клітинах A549.

Клітини A549 були піддані одноразовій дозі рентгенівського опромінення 2Гр або 4Гр за допомогою приладу X-RAD 160-225. Після опромінення клітини A549 (з одноразовою дозою опромінення), а також неопромінені клітини A549 були інкубовані ще 48 годин. Після цього було зібрано загальну РНК для високопродуктивного секвенування. Загальний дизайн експерименту: дослідження 3 умов по 3 репліки. Три групи: група 0Гр (контроль), група 2Гр та група 4Гр.

Відповідна стаття: doi: 10.1038/s41419-022-04561-x, PMID: 35190532.

Дані експерименту та метадані

Метадані зразків

kable(rse_colData_f |> 
  select(external_id, sra.study_title, sra.sample_title) |>  head(), format = "html", table.attr = "class='table table-striped'")
external_id sra.study_title sra.sample_title
SRR8371688 SRR8371688 RNA-sequencing in irradiated and normal A549 cells. 0GY-1
SRR8371689 SRR8371689 RNA-sequencing in irradiated and normal A549 cells. 0GY-2
SRR8371690 SRR8371690 RNA-sequencing in irradiated and normal A549 cells. 0GY-3
SRR8371691 SRR8371691 RNA-sequencing in irradiated and normal A549 cells. 2GY-1
SRR8371692 SRR8371692 RNA-sequencing in irradiated and normal A549 cells. 2GY-2
SRR8371693 SRR8371693 RNA-sequencing in irradiated and normal A549 cells. 2GY-3

Метадані генів

kable(head(rse_rowData_f), format = "html", table.attr = "class='table table-striped'")
source type bp_length phase gene_id gene_type gene_name level havana_gene tag
ENSG00000278704.1 ENSEMBL gene 2237 NA ENSG00000278704.1 protein_coding BX004987.1 3 NA NA
ENSG00000277400.1 ENSEMBL gene 2179 NA ENSG00000277400.1 protein_coding AC145212.2 3 NA NA
ENSG00000274847.1 ENSEMBL gene 1599 NA ENSG00000274847.1 protein_coding AC145212.1 3 NA NA
ENSG00000276256.1 ENSEMBL gene 2195 NA ENSG00000276256.1 protein_coding AC011043.1 3 NA NA
ENSG00000278198.1 ENSEMBL gene 1468 NA ENSG00000278198.1 protein_coding AC011043.2 3 NA NA
ENSG00000273496.1 ENSEMBL gene 1468 NA ENSG00000273496.1 protein_coding AC011841.1 3 NA NA

Матриця експресії

kable(rse_assay_f |> 
  arrange(desc(SRR8371688)) |> head(), format = "html", table.attr = "class='table table-striped'")
SRR8371688 SRR8371689 SRR8371690 SRR8371691 SRR8371692 SRR8371693 SRR8371694 SRR8371695 SRR8371696
ENSG00000198804.2 11438923 12624134 9876942 11286817 11500584 7659250 5559007 6564597 9430680
ENSG00000115414.18 6535629 7642890 5454817 6763706 5559684 4748784 3533072 3375268 4997642
ENSG00000198886.2 6206643 7357551 5660816 5944139 5889170 3652051 2882935 3346084 4630646
ENSG00000165092.12 5324460 5907729 4752635 6099076 6003428 4324184 2738698 2820161 4443411
ENSG00000198712.1 3950989 4529637 3600388 4185883 4218314 2704965 2025084 2605506 3313521
ENSG00000198786.2 3844190 4525495 3422576 3768348 3622434 2426284 1887319 2046555 3086360

Попередня візуалізація даних

Побудовано графік щільності для окремого зразка (SRR8371688) та гістограму частоти значень експресії

# Створення таблиці метаданих
colData_new <- data.frame(
  row.names = colnames(rse_filtered),
  condition = factor(c(rep("0Gy", 3), rep("2Gy", 3), rep("4Gy", 3)))
)

# Матриця counts
countData <- assay(rse_filtered)

# Побудова графіку щільності
p1 <- ggplot(countData, aes(x = SRR8371688 )) +
  geom_density(color = "darkblue", fill = "lightblue") +
  xlim(0, 10e4)

count_summary <- countData |>
  as.vector() |>
  data.frame(counts = _) |>
  group_by(counts) |>
  summarise(frequency = n(), .groups = 'drop') |>
  arrange(counts)


# Створення гістограми частоти значень
p2 <- ggplot(count_summary[1:200,], aes(x = counts, y = frequency, fill = frequency)) +
  geom_bar(stat = "identity", width = 0.7) +
  scale_y_log10(labels = scales::comma) +
  labs(title = "Histogram of Counts vs Frequency (Logarithmic Scale)",
       x = "Counts",
       y = "Frequency (log10)") +
  scale_fill_gradient(low = "#AA05BA", high = "#56B1F7") +
  theme_minimal()
#p1
#p2
counts_info <- p1 + p2

ggsave("./figs/counts_info.png", plot = counts_info, width = 16, height = 6, dpi = 300)
knitr::include_graphics("./figs/counts_info.png")

Аналіз диференційної експресії генів

DESeq

Для аналізу диференційної експресії генів було застосовано DESeq Wald

Коригування p.adj

Додатково було проведено коригування на хибний позитивний результат p.adj (FDR). Нижче наведено відповідний графік

res_wald_2Gy_vs_0Gy <- results(dds_wald, contrast = c("condition", "2Gy", "0Gy"))
res_wald_4Gy_vs_0Gy <- results(dds_wald, contrast = c("condition", "4Gy", "0Gy"))
res_wald_4Gy_vs_2Gy <- results(dds_wald, contrast = c("condition", "4Gy", "2Gy"))

hist_data_before <- bind_rows(
  data.frame(padj = res_wald_2Gy_vs_0Gy$padj, Comparison = "2Gy vs 0Gy"),
  data.frame(padj = res_wald_4Gy_vs_0Gy$padj, Comparison = "4Gy vs 0Gy"),
  data.frame(padj = res_wald_4Gy_vs_2Gy$padj, Comparison = "4Gy vs 2Gy")
)

res_wald_2Gy_vs_0Gy$padj <- p.adjust(res_wald_2Gy_vs_0Gy$pvalue, method = "fdr")
res_wald_4Gy_vs_0Gy$padj <- p.adjust(res_wald_4Gy_vs_0Gy$pvalue, method = "fdr")
res_wald_4Gy_vs_2Gy$padj <- p.adjust(res_wald_4Gy_vs_2Gy$pvalue, method = "fdr")

hist_data_after <- bind_rows(
  data.frame(padj = res_wald_2Gy_vs_0Gy$padj, Comparison = "2Gy vs 0Gy"),
  data.frame(padj = res_wald_4Gy_vs_0Gy$padj, Comparison = "4Gy vs 0Gy"),
  data.frame(padj = res_wald_4Gy_vs_2Gy$padj, Comparison = "4Gy vs 2Gy")
)

deseq_Wald_padj_b <- ggplot(hist_data_before, aes(x = padj, fill = Comparison, color = Comparison)) +
  geom_density(alpha = 0.2, adjust = 1.5) + # Криві густини
  scale_fill_manual(values = c("2Gy vs 0Gy" = "red", "4Gy vs 0Gy" = "blue", "4Gy vs 2Gy" = "forestgreen")) +
  scale_color_manual(values = c("2Gy vs 0Gy" = "red", "4Gy vs 0Gy" = "blue", "4Gy vs 2Gy" = "forestgreen")) +
  labs(title = "Розподіл коригованих p-значень", x = "p-значення", y = "Щільність") +
  geom_vline(xintercept = 0.05, linetype = "dashed", color = "black") +
  annotate("text", x = 0.05, y = Inf, label = "0.05", vjust = 1.5, hjust = -0.1, color = "black") +
  theme_minimal()

deseq_Wald_padj_a <- ggplot(hist_data_after, aes(x = padj, fill = Comparison, color = Comparison)) +
  geom_density(alpha = 0.2, adjust = 1.5) + # Криві густини
  scale_fill_manual(values = c("2Gy vs 0Gy" = "red", "4Gy vs 0Gy" = "blue", "4Gy vs 2Gy" = "forestgreen")) +
  scale_color_manual(values = c("2Gy vs 0Gy" = "red", "4Gy vs 0Gy" = "blue", "4Gy vs 2Gy" = "forestgreen")) +
  labs(title = "Розподіл коригованих p-значень", x = "p-значення", y = "Щільність") +
  geom_vline(xintercept = 0.05, linetype = "dashed", color = "black") +
  annotate("text", x = 0.05, y = Inf, label = "0.05", vjust = 1.5, hjust = -0.1, color = "black") +
  theme_minimal()

deseq_Wald_padj <- deseq_Wald_padj_b + deseq_Wald_padj_a
ggsave("./figs/deseq_Wald_padj.png", plot = deseq_Wald_padj, width = 16, height = 6, dpi = 300)
knitr::include_graphics("./figs/deseq_Wald_padj.png")

Volcano plot

Для візуалізації у зміні експресії генів для випадку опромінення рентгеном 2Гр (порівняння з контролнм) та 4Гр (порівняння з контролнм) побудовано volcano plot. Окремо проводилось порівняння експресії генів при опроміненні 2Гр та 4Гр.

# png("./figs/deseq_Wald.png", width = 2400, height = 800)
# par(mfrow = c(1, 3), mar = c(5, 5, 4, 2))

# plotMA(res_wald_2Gy_vs_0Gy, ylim=c(-30,30), main="MA-plot: 0Gy vs 2Gy")
# plotMA(res_wald_4Gy_vs_0Gy, ylim=c(-30,30), main="MA-plot: 0Gy vs 4Gy")
# plotMA(res_wald_4Gy_vs_2Gy, ylim=c(-30,30), main="MA-plot: 2Gy vs 4Gy")
# dev.off()
# knitr::include_graphics("./figs/deseq_Wald.png")

res_wald_df_2Gy_vs_0Gy <- as.data.frame(res_wald_2Gy_vs_0Gy)
res_wald_df_2Gy_vs_0Gy$significant <- res_wald_df_2Gy_vs_0Gy$padj < 0.05

res_wald_df_4Gy_vs_0Gy <- as.data.frame(res_wald_4Gy_vs_0Gy)
res_wald_df_4Gy_vs_0Gy$significant <- res_wald_df_4Gy_vs_0Gy$padj < 0.05

res_wald_df_4Gy_vs_2Gy <- as.data.frame(res_wald_4Gy_vs_2Gy)
res_wald_df_4Gy_vs_2Gy$significant <- res_wald_df_4Gy_vs_2Gy$padj < 0.05

res_wald_df_2Gy_vs_0Gy$comparison <- "2Gy vs 0Gy"
res_wald_df_4Gy_vs_0Gy$comparison <- "4Gy vs 0Gy"
res_wald_df_4Gy_vs_2Gy$comparison <- "4Gy vs 2Gy"

res_wald_df_2Gy_vs_0Gy <- as.data.frame(res_wald_df_2Gy_vs_0Gy)  |> 
  mutate(gene_id = rownames(res_wald_df_2Gy_vs_0Gy), comparison = "2Gy vs 0Gy")  |> 
  left_join(gene_annotations, by = "gene_id")

# Оновлення res_df_4Gy_vs_0Gy
res_wald_df_4Gy_vs_0Gy <- as.data.frame(res_wald_df_4Gy_vs_0Gy)  %>% 
  mutate(gene_id = rownames(.), comparison = "4Gy vs 0Gy")  |> 
  left_join(gene_annotations, by = "gene_id")

# Оновлення res_df_4Gy_vs_2Gy
res_wald_df_4Gy_vs_2Gy <- as.data.frame(res_wald_df_4Gy_vs_2Gy)  |> 
  mutate(gene_id = rownames(res_wald_df_4Gy_vs_2Gy), comparison = "4Gy vs 2Gy")  |> 
  left_join(gene_annotations, by = "gene_id")

# Побудова volcano plot
plot_2Gy_vs_0Gy_W <- EnhancedVolcano(res_wald_df_2Gy_vs_0Gy,
                                   lab = res_wald_df_2Gy_vs_0Gy$gene_name,
                                   x = 'log2FoldChange',
                                   y = 'padj',
                                   xlim = c(-5, 5),
                                   ylim = c(0, 10),
                                   pCutoff = 0.05,
                                   FCcutoff = 1.5,
                                   pointSize = 5,
                                   labSize = 6,
                                   col = c('grey30', 'forestgreen', 'royalblue', 'red2'),
                                   title = '2Gy vs 0Gy',
                                   legendPosition = 'top')

plot_4Gy_vs_0Gy_W <- EnhancedVolcano(res_wald_df_4Gy_vs_0Gy,
                                   lab = res_wald_df_4Gy_vs_0Gy$gene_name,
                                   x = 'log2FoldChange',
                                   y = 'padj',
                                   xlim = c(-5, 5),
                                   ylim = c(0, 10),
                                   pCutoff = 0.05,
                                   FCcutoff = 1.5,
                                   pointSize = 5,
                                   labSize = 6,
                                   col = c('grey30', 'forestgreen', 'royalblue', 'red2'),
                                   title = '4Gy vs 0Gy',
                                   legendPosition = 'top')

plot_4Gy_vs_2Gy_W <- EnhancedVolcano(res_wald_df_4Gy_vs_2Gy,
                                   lab = res_wald_df_4Gy_vs_2Gy$gene_name,
                                   x = 'log2FoldChange',
                                   y = 'padj',
                                   xlim = c(-5, 5),
                                   ylim = c(0, 10),
                                   pCutoff = 0.05,
                                   FCcutoff = 1.5,
                                   pointSize = 5,
                                   labSize = 6,
                                   col = c('grey30', 'forestgreen', 'royalblue', 'red2'),
                                   title = '4Gy vs 2Gy',
                                   legendPosition = 'top')
volc_wald_plot <- plot_2Gy_vs_0Gy_W + plot_4Gy_vs_0Gy_W + plot_4Gy_vs_2Gy_W + plot_layout(nrow = 1)

#print(volc_wald_plot)
ggsave("./figs/volc_wald.png", plot = volc_wald_plot, width = 20, height = 6, dpi = 300)

2Gy vs 0Gy

ggsave("./figs/volc_plot_2Gy_vs_0Gy.png", plot = plot_2Gy_vs_0Gy_W, width = 20, height = 20, dpi = 300)
knitr::include_graphics("./figs/volc_plot_2Gy_vs_0Gy.png")

4Gy vs 0Gy

ggsave("./figs/volc_plot_4Gy_vs_0Gy.png", plot = plot_4Gy_vs_0Gy_W, width = 20, height = 20, dpi = 300)
knitr::include_graphics("./figs/volc_plot_4Gy_vs_0Gy.png")

4Gy vs 2Gy

ggsave("./figs/volc_plot_4Gy_vs_2Gy.png", plot = plot_4Gy_vs_2Gy_W, width = 20, height = 20, dpi = 300)
knitr::include_graphics("./figs/volc_plot_4Gy_vs_2Gy.png")

Upset

Для візуалізації перекриття між списками генів, які є статистично значущими (з коригованим p-значенням < 0.05) та мають великий зсув логарифмічних змін (log2FoldChange > 1.5 для підвищеної експресії та log2FoldChange < -1.5 для зниженої експресії) у різних порівняннях, використали upset plot.

res_2Gy_0Gy_up <- res_wald_df_2Gy_vs_0Gy |> as.data.frame() |> filter(pvalue < 0.05 & log2FoldChange > 1.5)
res_2Gy_0Gy_down <- res_wald_df_2Gy_vs_0Gy  |> as.data.frame() |> filter(pvalue < 0.05 & log2FoldChange < -1.5)
res_4Gy_0Gy_up <- res_wald_df_4Gy_vs_0Gy  |> as.data.frame() |> filter(pvalue < 0.05 & log2FoldChange > 1.5)
res_4Gy_0Gy_down <- res_wald_df_4Gy_vs_0Gy  |> as.data.frame() |> filter(pvalue < 0.05 & log2FoldChange < -1.5)
res_4Gy_2Gy_up <- res_wald_df_4Gy_vs_2Gy |> as.data.frame() |> filter(pvalue < 0.05 & log2FoldChange > 1.5)
res_4Gy_2Gy_down <- res_wald_df_4Gy_vs_2Gy  |> as.data.frame() |> filter(pvalue < 0.05 & log2FoldChange < -1.5)

venn_list <- list("2Gy vs 0Gy(up)" = res_2Gy_0Gy_up$gene_name,
                  "2Gy vs 0Gy(down)" = res_2Gy_0Gy_down$gene_name,
                  "4Gy vs 0Gy(up)" = res_4Gy_0Gy_up$gene_name,
                  "4Gy vs 0Gy(down)" = res_4Gy_0Gy_down$gene_name,
                  "4Gy vs 2Gy(up)" = res_4Gy_2Gy_up$gene_name,
                  "4Gy vs 2Gy(down)" = res_4Gy_2Gy_down$gene_name)

ven <- venndetail(venn_list)

png("./figs/upset_diagram.png", width = 1600, height = 1600, res = 300)
plot(ven, type = "upset")
dev.off()

Аналіз головних компонент (PCA)

Наступним кроком було проведено PCA аналіз. Scree plot Червона лінія відображає кумулятивну суму варіацій, яка показує, що 80% варіації даних припадає на п’ять головних компонент. При цьому PC1 пояснює близько 50% варіації. Для того, щоб визначити кількість головних компонент, які варто зберегти для аналізу були використані метод Горна та elbow, які зображені на графіку штриховою лінією. Biplot комбінує інформацію про зразки і змінні в одному графіку. По осях зображені дві головні компоненти: PC1 і PC2, кожна з яких пояснює відсоток варіації в даних: PC1 — 44.98%, PC2 — 12.94%. Точки відповідають різним спостереженням за умов випромінювання X-хвилями 0Гр, 2Гр, 4Гр, які приблизно кластеризуються відповідно умовам впливу. Щоб перевірити чи є “правильна” кластеризація експериментів, був проведений кластерний аналіз методом k-means з 5-ма головними компонентами, адже саме стільки компонент пояснюють більшу варіацію даних. Loadings plot показуює, наскільки кожен ген впливає на кожну головну компоненту. Це значення визначає, наскільки цей ген сприяє поясненню варіації в даних.

vst_data <- assay(vst(dds_wald))

pca_result <- PCAtools::pca(vst_data, metadata = colData_new, removeVar = 0.1)

horn <- parallelPCA(vst_data)
elbow <- findElbowPoint(pca_result$variance)
var_explained <- which(cumsum(pca_result$variance)> 80)[1]

pscree <- PCAtools::screeplot(pca_result,
    axisLabSize = 18,
    titleLabSize = 22,
    vline = c(horn$n, elbow, var_explained), hline = 85,
    returnPlot = FALSE) +
    geom_label(aes(x = horn$n, y = 60,
      label = 'Horn\'s', vjust = -1, size = 8)) +
    geom_label(aes(x = elbow, y = 60,
      label = 'Elbow method', vjust = -3, size = 8)) +
    geom_label(aes(x = var_explained , y = 85,
      label = '80% explained variation', vjust = -1, size = 8)) 


ppairs <- PCAtools::pairsplot(pca_result,
    hline = 0, vline = 0,
    pointSize = 2,
    gridlines.major = FALSE, gridlines.minor = FALSE,
    colby = 'condition',
    plotaxes = TRUE,
    margingaps = unit(c(0.01, 0.01, 0.01, 0.01), 'cm'),
    returnPlot = FALSE)


# Створюємо нові назви для рядків, зіставляючи gene_id з gene_name
new_row_names <- gene_annotations$gene_name[match(rownames(pca_result$loadings), gene_annotations$gene_id)]

# Якщо є дублікат імен, додамо індекси для унікальності
new_row_names <- make.unique(new_row_names)

# Перезаписуємо назви рядків у pca_result$loadings
rownames(pca_result$loadings) <- new_row_names


pbiplot <- PCAtools::biplot(pca_result,
    axisLabSize = 18,
    # loadings parameters
    showLoadings = TRUE,
    sizeLoadingsNames = 4,
    colLoadingsNames = 'red4',
    #points
    lab = NULL,
    colby = 'condition',
    colLegendTitle = 'Condition',
    colkey = c('0Gy' = '#228B22', '2Gy' = '#7B68EE', '4Gy' = '#FF8C00'),
    pointSize = 5,
    #plotview
    hline = 0,
    vline = c(-5, 0, 5),
    vlineType = c('dotdash', 'solid', 'dashed'),
    gridlines.major = FALSE, gridlines.minor = FALSE,
    legendPosition = 'right', 
    legendLabSize = 16,
    legendIconSize = 8.0,,
    subtitle = 'PC1 versus PC2',
    caption = '5 PCs ≈ 80%',
    returnPlot = FALSE)

ploadings <- PCAtools::plotloadings(pca_result,
    axisLabSize = 18,
    subtitle = 'PC1, PC2, PC3, PC4, PC5',
    labSize = 3,
    returnPlot = FALSE)

pca_df <- as.data.frame(pca_result$rotated)

# Додаємо інформацію про condition та інші метадані
pca_df$condition <- colData_new$condition

# Використовуємо більше компонент для кластеризації, наприклад, перші 5
set.seed(1234)
km <- kmeans(x = pca_df[, 1:5], centers = 3, nstart = 25)  # Використовуємо 5 компонент і 25 запусків
# Додаємо результати кластеризації в датафрейм
pca_df$cluster <- as.factor(km$cluster)

# Візуалізація
pkmean <- autoplot(km, data = pca_df, frame = TRUE, frame.type = 'norm', size = 4, alpha = 0.6) + 
  aes(shape = condition, color = cluster) + 
  scale_shape_manual(values = c('0Gy' = 16, '2Gy' = 17, '4Gy' = 18)) + 
  labs(subtitle = 'with 5 PCs',
       color = 'cluster', shape = 'Condition') +
  coord_equal() + 
  theme_minimal(base_size = 18)+
  theme(legend.position = 'right', # Позиція легенди
        panel.border = element_rect(color = 'black', fill = NA, size = 1))  # Додання рамки


# Додаємо заголовки до кожного графіка
pscree <- pscree + ggtitle("A Scree plot") + theme(plot.title = element_text(face = "bold", hjust = -0.01, size = 18))
ploadings <- ploadings + ggtitle("B Loading plot") + theme(plot.title = element_text(face = "bold", hjust = -0.01, size = 18))
pbiplot <- pbiplot + ggtitle("C Bi-plot") + theme(plot.title = element_text(face = "bold", hjust = -0.01, size = 18))
pkmean <- pkmean + ggtitle("D K-means Clusters") + theme(plot.title = element_text(face = "bold", hjust = -0.01, size = 18))
ppairs <- ppairs + ggtitle("E Pairs plot") + theme(plot.title = element_text(face = "bold", hjust = 0.5, size = 18))

# Комбінування графіків у сітку
combined_plot <- (pscree | ploadings) / (pbiplot | pkmean) / ppairs


ggsave("./figs/combined_plot.png", combined_plot, width = 20, height = 30)

knitr::include_graphics("./figs/combined_plot.png")

Функціональний аналіз

Після проведення аналізу диференційної експресії було визначено роль цих генів в клітинах. Для цього було проведено аналіз збагачення за допомогою Enrichr

Функціональний аналіз Enrichr

Узагальнене

knitr::include_graphics("./figs/enrichr_results_all_up.png")

2Gy 0Gy

knitr::include_graphics("./figs/enrichr_results_2Gy_0Gy_up.png")

4Gy 0Gy

knitr::include_graphics("./figs/enrichr_results_4Gy_0Gy_up.png")

4Gy 2Gy

knitr::include_graphics("./figs/enrichr_results_4Gy_2Gy_up.png")

Аналіз метаболічних шляхів за допомогою Pathview

Для того, щоб явно побачити роль генів з підвищеною та зі зниженою експресією у клітинах лінії А549 було проведено аналіз метаболічних шляхів циклу клітини, проліферації, апоптозу.

logFC <- path_4Gy_vs_0Gy$log2FoldChange
names(logFC) <- path_4Gy_vs_0Gy$entrezgene_id


pathview(gene.data = logFC, 
         pathway.id = "hsa05207", 
         species = "hsa", 
         #kegg.native = T
         ut.suffix = "Non-small cell lung cancer (специфічно для A549)")

pathview(gene.data = logFC, 
         pathway.id = "hsa04110", 
         species = "hsa", 
         #kegg.native = T
         ut.suffix = "Cell cycle")

pathview(gene.data = logFC, 
         pathway.id = "hsa04115", 
         species = "hsa", 
         #kegg.native = T
         ut.suffix = "p53 signaling pathway")

pathview(gene.data = logFC, 
         pathway.id = "hsa04725", 
         species = "hsa", 
         #kegg.native = T
         ut.suffix = "NAChR")

pathview(gene.data = logFC, 
         pathway.id = "hsa04210", 
         species = "hsa", 
         #kegg.native = T
         ut.suffix = "Apoptosis")

Канцерогенез

knitr::include_graphics("hsa05207.pathview.png")

Клітинний цикл

knitr::include_graphics("hsa04110.pathview.png")

р53 сигнальний шлях

knitr::include_graphics("hsa04115.pathview.png")

Апоптоз

knitr::include_graphics("hsa04210.pathview.png")

Холінергічний синапс

knitr::include_graphics("hsa04725.pathview.png")

## {-}